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Iterative LQG Controller Design 
Through Closed-Loop Identification 

Min-Hung Hsiao , 1 Jen-Kuang Huang , 2 and 
David E. Cox 3 


This paper presents an iterative LQG controller design ap- 
proach for a linear stochastic system with an uncertain open- 
loop model and unknown noise statistics. This approach consists 
of closed-loop identification and controller redesign cycles. In 
each cycle, the closed-loop identification method is used to 
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identify an open-loop model and a steady-state Kalman filter 
gain from closed-loop input/output test data obtained by using 
a feedback LQG controller designed from the previous cycle. 
Then the identified open-loop model is used to redesign the 
state feedback. The state feedback and the identified Kalman 
filter gain are used to form an updated LQG controller for the 
next cycle. This iterative process continues until the updated 
controller converges. The proposed controller design is demon- 
strated by numerical simulations and experiments on a highly 
unstable large-gap magnetic suspension system. 


1 Introduction 

Classical Linear Quadratic Gaussian (LQG) controllers are 
designed by solving two separate, but dual problems: the Linear 
Quadratic Regulator (LQR) design and Kalman filter (i.e., opti- 
mal state estimator) design. The performance of the controllers 
relies on an accurate open-loop model for the LQR and an 
accurate estimate of the measurement and process noise statis- 
tics for the Kalman filter. It is difficult to obtain an accurate 
model through analysis for some systems, and an accurate esti- 
mate of the noise statistics through testing for most systems. 
Furthermore, the noise statistics may be related to the controller 
if part of the measurement and process noise are generated by 
the sensor and actuator amplifiers, respectively. To overcome 
these problems, we present an iterative LQG controller design 
approach for a linear stochastic system with an uncertain open- 
loop model and unknown noise statistics. This approach consists 
of closed-loop identification and controller redesign cycles. The 
closed-loop identification method can simultaneously identify 
the open-loop model and the Kalman filter gain under the 
closed-loop operation with a known dynamic controller. Then 
the identified open-loop model is used for the LQR design. The 
LQR and the identified Kalman filter gain are used to form the 
updated LQG controller for the next closed-loop identification. 
The process continues until the updated LQG controller con- 
verges. 

For system identification, several methods (Chen et al., 
1992a, 1992b, 1993; Phan et aL, 1991; Juang et al., 1993) have 
been introduced recently to identify the state-space model of a 
linear system and the Kalman filter. Typically the system is 
under open-loop excitation with an uncorrelated white noise 
input. For an unstable system, the input/output data are not 
available while it is under an open-loop operation. To directly 
use these methods, we have to design a controller and an input 
signal for the closed-loop system so that the input signal to the 
open-loop system is almost white. Unfortunately, this is very 
difficult. On the other hand, some identification methods (Phan 
et al., 1992; Liu and Skelton, 1990) have been proposed recently 
for identifying a system under closed-loop operation. However, 
they have several shortcomings. First, the Kalman filter gain 
cannot be simultaneously identified because they are applied 
only for deterministic systems. In Phan et al. (1992), no re- 
cursive form was derived for computing the open-loop system 
Markov parameters, and in Liu and Skelton (1990), the ap- 
proach is based on system pulse response. In this paper, a re- 
cursive form for computing the open-loop system and Kalman 
filter Markov parameters is derived for stochastic systems with 
random excitation. 

For a system under closed-loop operation, a novel approach 
for identifying the open-loop model and Kalman filter gain is 
presented. First, we derive the relation between closed-loop 
state-space and AutoRegressive with eXogeneous ( ARX) mod- 
els for stochastic systems. From the derivation, it can be seen 
that a state-space model can be represented by an ARX model 
if the order of the ARX model is chosen large enough. Since 
the relation between the input/output data and the system pa- 
rameters of an ARX model is linear, a linear programming 
approach like least-square methods, can be used for the ARX 
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model parameter estimation. Second, we derive the algorithm 
to compute the open-loop system and Kalman filter Markov 
parameters from the estimated ARX model parameters. In this 
step, we first compute the closed-loop system and Kalman filter 
Markov parameters from the estimated ARX model parameters. 
Then the open-loop system and Kalman filter Markov parame- 
ters are computed from the closed-loop system and Kalman 
filter Markov parameters and the known controller Markov pa- 
rameters. Third, the state-space model for the open-loop system 
is realized from the open-loop Markov parameters through the 
singular value decomposition method (Chen et al., 1984; Juang 
and Pappa, 1985). Finally, the Kalman filter for the open-loop 
system can be estimated from the realized state-space model 
and the open-loop Kalman filter Markov parameters through a 
least-square approach. 

With this closed-loop identification, an iterative LQG control- 
ler design can be performed. Since the Kalman filter used in 
this LQG controller is obtained directly from the closed-loop 
identification, it automatically takes into account the effect of 
the controller on the noise statistics. The LQR tends to reject 
the process noise and the Kalman filter tends to filter out the 
measurement noise. Therefore, the closed-loop identification 
can improve the LQG design and an updated LQG controller 
can enhance the closed-loop identification in the next cycle. 
After a certain number of iterations, the LQG controller will 
converge. 

A similar approach is presented by Liu and Skelton ( 1990). 
As compared to that approach, this paper has the following 
contributions. First, the proposed method is developed under 
the stochastic framework rather than a deterministic one. Sec- 
ond, the Kalman filter gain is also identified so that it can be 
used for state estimation directly. Third, random excitation 
rather than pulse response is used for the closed-loop identifica- 
tion. Finally, since the Kalman filter gain is identified, LQR 
state feedback is used rather than output feedback. Numerical 
and experimental examples are provided to illustrate and vali- 
date this controller design. 

2 Closed-Loop State-Space and ARX Models Rela- 
tionship 

Since the relation between the input/output data and the 
model parameters of a state-space model is nonlinear, parameter 
estimation of a state-space model from input/output data is 
a nonlinear programming problem. Nonlinear programming is 
difficult to solve in general and involves complex iterative nu- 
merical methods. The convergence and uniqueness of the solu- 
tion are also not guaranteed. Unlike a state-space model, the 
ARX model has a linear relationship between its model parame- 
ters and input/output data. Therefore, linear programming can 
be used for identifying the ARX model. After obtaining the 
ARX model, a state-space model can be developed based on the 
relation between these two models. In this section, the relation 
between a closed-loop state-space model and an ARX'Thodel is 
derived by using z-transforms. 

A finite-dimensional, linear, discrete-time, time-invariant 
system can be modeled as: 

x* + , = Ax k + Bu k + w* (1) 

y k = Cx k + v k (2) 

where * € R nX 1 , « € R sX 1 , y G R mX 1 are state, input, and output 
vectors, respectively; w* is the process noise, v* the measurement 
noise; [A, B, C] are the state-space parameters. Sequences w k 
and v* are assumed gaussian, white, zero-mean, and stationary 
with covariance matrices W and V, respectively. One can derive 
a steady-state filter innovation model (Hay kin, 1991): 

x k+1 = Ax k + Bu k + AKe k (3) 

y k = Cx k + e*. (4) 
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where x k is the a priori estimated state, K is the steady-state 
Kalman filter gain, and e* is the residual after filtering. The 
existence of K is guaranteed if the system is detectable and (A , 
W U1 ) is stabilizable (Goodwin and Sin, 1984). The advantage 
of using the filter innovation model in the closed-loop identifi- 
cation is that one can directly identify the Kalman filter gain 
without estimating the covariance matrices of both process and 
measurement noise which usually are difficult to be obtained 
from test data (Chen and Huang, 1994). 

Similarly, any kind of dynamic output feedback controller 
can be modeled as: 

P*+i = A d p k + B d y k . (5) 

= C d p k + D d y k + r kt (6) 

where A d ,B d ,C d , and D d are the system matrices of the dynamic 
output feedback controller, p k the controller state and r k is the 
open-loop input to the closed-loop system. Combining (3) to 
(6), the augmented closed-loop system dynamics becomes 
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b g ] T , and / is the number of data points. The integer l has to 
be large enough so that the 4> matrix has more rows than col- 
umns. The batch least-square solution is 

$ = ($ T Qy l $ T Z (17) 

Therefore, solving for an ARX model simply involves solving 
a linear programming problem involving an over determined 
set of equations. 

3 Markov Parameters and State-Space Realization 


77* +J = A c r} k + B c r k + A e K c e k 


where 


yk = C c r) k + e*, 



(7) 

( 8 ) 


(9) 


and j) k = [*»]. It is noted that K c can be considered as the 
Kalman filter gain for the closed-loop system and the existence 
of the steady-state K c is guaranteed when the closed-loop system 
matrix A c is nonsingular. Substituting (8) into (7) yields 

77* +i = Ai) k + B c r k + A c K c y k , (10) 

where A - A c — A c K c C c and is guaranteed to be asymptotically 
stable because the steady-state Kalman filter gain K c exists. The 
z-transform of (10) and (8) yields 


In the previous section, an ARX model, which represents a 
closed-loop system, is identified from the input /output data 
through the least-square method. With the known controller 
dynamics, the estimated ARX model can be transformed to an 
open-loop state-space model by the following steps. First, the 
closed-loop system and Kalman filter Markov parameters are 
calculated from the estimated coefficient matrices of the ARX 
model. Second, the open-loop system and Kalman filter Markov 
parameters are derived from the closed-loop system Markov 
parameters, the closed-loop Kalman filter Markov parameters, 
and the known controller Markov parameters. Third, the open- 
loop state-space model is realized by using singular-value de- 
composition for a Hankel matrix formed by the open-loop sys- 
tem Markov parameters. Finally, an open-loop Kalman filter 
gain is calculated from the realized state-space model and the 
open-loop Kalman filter Markov parameters through least- 
squares. 

The z-transform of the open-loop state-space model (3 ) yields 

x (z) = (z - A)- l (Bu(z) + AKe(z)). (18) 

Substituting ( 18) to the z-transform of the output Eq. (4), one 
has 


Viz) = (Z - A)~ l (A e K c y(z) + B c r(z)) (11) 
y(z) = C c tj(z) + e(z). (12) 

Substituting (11) into (12), one has 

yiz) = C f (z - A)~ l (A e K c y(z) + B c r ( Z )) + e(z). (13) 

The inverse z-transform of (13) with (z - A) -1 = 2 A'~ l z~‘ 

f- 1 

yields 

OD 00 

y* = 2 C c S l - x A e K e y k . l + I C c A‘- l B c r k ., + e k . (14) 

i- 1 >-l 

Since A is asymptotically stable, A‘ « 0 if i > q for a sufficient 
large number q. Thus (14) becomes 

? <7 

y* « 2 a t y k -i + 2 V*-/ + e* (15) 

i-i »— i 

where 

a, = C c A‘~ l A c Kc, bi = C c X f - l B c . (16) 


y(z) = C(z - A)~'(Bu(z) + AKe(z)) + e(z) 

= 2 Y(k)z~ k u(z ) + 2 N(k)z~ k e(z), (19) 

*-i *«o 

where Y (k) = CA k ~'B are the open-loop system Markov pa- 
rameters; N(k) = CA k ~ l AK, for k = 1 , . . . , open-loop 
Kalman filter Markov parameters, and N{ 0) = I which is an 
identity matrix. Similarly, for the dynamic output feedback con- 
troller (5) and (6) and the closed-loop state-space model (7) 
and (8), one can derive ^ 

uiz) = 2 Y d (k)z~ k y(z) + r(z) (20) 

*-0 

00 go 

yiz) = 2 Y c (k)z~ k r(z) + 2 N c ik)z- k e(z). (21) 

*-t km 0 

where Y d ( 0) = D d , and Y d (k) = C d A d ~ l B d , for k = 1 ®, 

are the controller Markov parameters; Y c (k) = C c A k ~ l B c the 
closed-loop system Markov parameters; and N c ( 0) = 7, N c (k) 
— C c A k 1 A C K C , for k — 1, . . . , t», the closed-loop Kalman filter 
Markov parameters. 


The model described by (15) is the ARX model which di- 
rectly represents the relationship between the input and output 
of the closed-loop system. The coefficient matrices a, and b, 
can be estimated through least-square methods from random 
excitation input r k and the corresponding output y k . From (15) 


Closed-Loop System and Kalman Filter Markov Parame- 
ters. The z-transform of the ARX model (15) yields 
q q 

(7-2 aiZ'^yiz) = 2 biZ~Hz) + e(z). (22) 

i-l i— | 
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Applying long division to (22), one has 

y(z) = (bi z~ l + (b 2 + a\b x )z~ 2 + ( b 3 + a x {b 7 + a x b x ) 

+ a 2 b x )z ~ 3 + . . . )r(z) + (/ + a x z~' + (a,a, + a 2 )z~ 2 

+ (a 1 (fl,a 1 + a 2 ) + + a 3 )z~ 3 + . . . )c(z). 

After comparing with (21 ) , the closed-loop system and Kalman 
filter Markov parameters can be recursively calculated from the 
estimated coefficient matrices of the ARX model, 

k 

Y c {k) = b k + X aMk - i) (23) 

I- I 

* 

K(k) = I a,N e (k - i). (24) 

It is noted that y f (0) = 0, N c (0) = /, and a, = Z>, = 0, when 
i> q. One may obtain (23) and (24) from ( 16) and the defini- 
tion of the Markov parameters (Phan et al., 1991; Juang et al., 
1993). However, the derivation is much more complex. 

Open-Loop System and Kalman Filter Markov Parame- 
ters. Next, the open-loop system and Kalman filter Markov 
parameters can be derived from the closed-loop system Markov 
parameters, the closed-loop Kalman filter Markov parameters, 
and the known controller Markov parameters. Substituting (20) 
into ( 19) yields 

00 tlk 

y(z) = (X y(*)z-*)(I Y d (k)z- k y(z )) 

*- I k-0 

® 00 9 

+ X Y(k)z~ k r(z) + X N{k)z~ k t(z)= X a*z"*y(z) 

*-i *-o *-i 

+ X Y(k)z~ k r(z) + X N(k)z~ k e(z), (25) 

oo 

where a k = 1 Y (i)Y d (k - /). Rearranging (25), one has 

f-1 

oo 

(/ - X Qf*2“*)y(z) 

<E 00 

= X Y(k)z~ k r(z) + X N(k)z~ k e(z). (26) 

*-i *-o 


Equations (29) and (30) show that one can recursively calculate 
the open-loop system and Kalman filter Markov parameters 
from the closed-loop system Markov parameters in (23), the 
closed-loop Kalman filter Markov parameters in (24), and the 
known controller Markov parameters Y d (k) = C d A k d ~'B d . It is 
noted that K c (0) = 0 and N c ( 0) = I. One can easily verify (29) 
and (30) from (9), and also from the definition of the Markov 
parameters. 


State-Space Realization. The open-loop state-space model 
can be realized from the open-loop system Markov parameters 
through the Singular Value Decomposition (SVD) method 
(Chen, 1984; Juang and Pappa, 1985). The first step is to form 
a Hankel matrix from the open-loop system Markov parameters. 


HU) = 


YU) 

YU+D 


YU+ 1) 
Y(j + 2) 


Y(j + 0) 
YU + 0+1) 


Lyu + y) yu + 7+D 


yu + y + 0). 


(31) 

where Y(j) is the y-th Markov parameter. For a noise free 
system, if the arbitrary integers 0 & n, and y n (the order 
of the system), the Hankel matrix H(j) is of rank n. From the 
measurement Hankel matrix, the realization uses the SVD of 
H(l), H(l) = ULV T , to identify a n-th order discrete state- 
space model as 


A = l.: in UlH{2)VXn Xn , 

B = V n ‘ 2 V T n E„C = E T m U n iy 2 (32) 

where matrix X„ is the upper left hand n X n partition of Z 
containing the n largest singular values along the diagonal. Ma- 
trices U„ and V„ are obtained from U and V by retaining only 
the n columns of singular vectors associated with the n singular 
values. Matrix E m is a matrix of appropriate dimension having 
m columns, all zero except that the top m x m partition is an 
identity matrix. E, is defined similarly. 


Open-Loop Kalman Filter Gain. Once the open-loop A 
and C are obtained, one can easily calculate the open-loop Kal- 
man filter gain from the open-loop Kalman filter Markov param- 
eters A^(Ar) = CA k K in a least-square sense as follows 


Similarly, one can apply long division to (26), and then com- 
pare it with (21), to describe the closed-loop system Markov 
parameters recursively in terms of the open-loop system and 
the controller Markov parameters. 

r,U) = YU) + X a,Y,U - k) 

km 1 

j * 

= Y(j)+HY(i)Y d (k-i)Y c (j~k) (27) 

*»i i-i 

And the closed-loop Kalman filter Markov parameters can be 
recursively expressed in terms of the open-loop system Markov 
parameters, the open-loop Kalman filter Markov parameters, 
and the controller Markov parameters as follows: 

A r Aj) = Nij) + X ct k N c (j ~ k) 

km 1 

J k 

= N(j) +11 Y (i)Y d (k - i)N c (j - k). (28) 

*«] i— I 

Rearranging (27) and (28), one has 

j * ' 

YU) = YJj) - X X Y(i)YJk - i)Y, (j - k) (29) 

i-l i-J 

j * 

N(j) = N c (j) -11 Y(i)Y d (k - i)N c (j - k). (30) 

km 1 |-| 
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mN(k)~ 


_ CA k _ 


The integer k has to be large enough so that the matrix O has 
more rows than columns. The identified Kalman filter gain can 
be used directly for state estimation. 

4 Iterative LQG Controller Design __ 

Classical LQG controllers are designed by solving two sepa- 
rate, but dual problems: the LQR design and Kalman filter 
design. Here, the Kalman filter gain can be simultaneously ob- 
tained with the open-loop state-space model through the closed- 
loop identification. Only the LQR design based on the identified 
open-loop model needs to be solved. The performance index 
for the LQR problem is defined as 

oc 00 

P.I. = 1 ylQyk + ulRu k = X x\C T QCx k + u T k Ru k (34) 

*-l km\ 

where weighting matrices Q and R are design parameters. We 
can summarize the iterative LQG controller design as follows: 

1. Use the a priori open-loop model and arbitrary covari- 
ance matrices of the measurement and process noise to design 
the LQR and Kalman filter. Then, calculate the controller Mar- 
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Fig- 1 Large-angle magnetic suspension test facility (LAMSTF) 
configuration 


kov parameters. The weighting matrices Q and R for the LQR 
chosen here will remain the same in the following iterations. 

2. Apply random excitation input to the closed-loop system 
and record the closed-loop input/ output data. 

3. Estimate the coefficient matrices of the closed-loop ARX 
model by using (17). 

4. Calculate the closed-loop system and Kalman filter Mar- 
kov parameters by using (23) and (24). 

5. Calculate the open-loop system and Kalman filter Mar- 
kov parameters by using (29) and (30). 

6^ Realize the open-loop state-space system matrices [A, 
B, C] by using (31) and (32). 

7. Estimate the open-loop Kalman filter gain K by using 
(33). 

8. Obtain the LQR feedback gain F by solving the corre- 
sponding Riccati equation based on the identified open-loop 
model. 

9. Form the updated LQG controller in (5) and (6) by 
using A d = A - BF ~ AKC , B d = AK, C d = —F, and D d = 0. 

1 0. Calculate the updated controller Markov parameters and 
check the convergence of the controller by 

n 

6 = I IllrfWupdu- - ^(fcWousIh- (35) 

*-o 

If <5 is greater than a desired value, go back to step 2, otherwise 
stop. 

5 Numerical and Experimental Example 

The proposed iterative LQG controller design has been ap- 
plied to control design of the Large-Angle Magnetic Suspension 
Test Facility (LAMSTF) (Groom and Britcher, 1992; Groom 
and Schaffner, 1990) developed in NASA Langley Research 
Center ( see Fig. 1 ) . The LAMSTF is a laboratory-scale research 
project to demonstrate the magnetic suspension of objects over 
wide ranges of attitudes. This system represents a scaled model 
of a planned Large-Gap Magnetic Suspension System. The 
LAMSTF system consists of a planar array of five copper elec- 
tromagnets which actively suspend a small cylinder with a per- 
manent magnet core. The cylinder is a rigid body and has six 
independent degrees of freedom, namely, three displacements 
(x, y, and z) and rotations (pitch, yaw, and roll). Currents in 
the electromagnets generate a magnetic field which produces a 
net force and torque on the suspended cylinder. The roll of the 
cylinder is uncontrollable, and is assumed to be motionless. 
Five pairs of the LEDs and light receivers are used to indirectly 
sense the pitch and yaw angles, and three displacements of the 
cylinder’s centroid. Therefore, the control inputs to the system 
consist of five currents sent into five electromagnets and the 
system outputs are five voltage signals measured from five 
photo detectors. The forces on the cylinder are, in general, non- 
linear functions of space and current. Therefore, only the linear 
time-invariant perturbed motion about an equilibrium state is 



considered. Because it is difficult to accurately model the mag- 
netic field and its gradients, the analytical model contains some 
modeling errors. Therefore, the performance of the LQG con- 
troller based on the analytical model alone is unsatisfactory. 

The system matrices of the analytical model are shown in 
the appendix. The eigenvalues of the system matrix indicate 
that the LAMSTF system includes highly unstable real poles 
(about 10 Hz) and low-frequency oscillatory modes (about 0.16 
Hz). For both numerical simulation and experiment, the sam- 
pling rate is 250 Hz. The performance index used for the LQR 
design is also shown in the Appendix. The step command for 
all simulations and experiments is 0.02 radian for pitch and 
yaw, and 0.2 mm for x, y, and z. 

In the numerical simulation, the analytical model is used as 
the true model. In each iteration, the ratios of the process and 
measurement noise to the corresponding signal are 2 and 1 
percent, respectively. To simulate modeling error and unknown 
noise statistics, the initial LQG controller is designed by using 
an altered model of which each parameter is 5 percent greater 
than the corresponding parameter of the analytical model and 
guessed covariance matrices of noise W = 10/| OX io and V = 
/ 3 x 5 • The simulated step response with this initial controller for 
the pitch, yaw, x, y, and z is shown in Fig. 2. It is clear that 
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the result is very poor. After performing the first iteration of 
the proposed iterative LQG controller design, the step response 
shown in Fig. 3 is greatly improved. The performance is further 
improved slightly in the following iterations. Figure 4 shows 
how the controller converges by comparing the (1, 1) element 
of the controller Markov parameters. 

For a noise free system, the exact open-loop model can be 
obtained after the first closed-loop identification and no further 
iteration is required. In this case, the identified Kalman filter 
gain becomes the dead-beat observer gain (Phan et al., 1991; 
Juang et al., 1993). For a noise corrupted system, iterations are 
required to update the open-loop model and the Kalman filter 
gain until the iterative LQG controller converges. Although 
the numerical simulations show that the iterative controller can 
converge quickly, the required conditions to guarantee the con- 
vergence need further study. 

In the experiments, the analytical model and guessed covari- 
ance matrices of noise W = 10/, O xio and V — l 5x5 are used to 
design the initial LQG controller. The experimental step re- 
sponse with this initial controller is also very poor. The experi- 
mental step responses for the first three iterations are compared 
in Fig. 5 to demonstrate how the step response is improved with 
iteration. In each iteration, the open-loop system model and the 
Kalman filter gain are updated through the closed-loop identifi- 
cation from experimental data. The experimental step response 
improves with each iteration, similar to the simulated cases. 
The experimental steady-state errors, however, do not go to 
zero in each case. This is due to drift in the sensor zero between 
experiments. The system’s dynamics have been found to be 
insensitive to these small changes in the operating point. The 
results show that the proposed iterative LQG controller design 
is very effective for controlling this highly unstable magnetic 
suspension system. 




Isiitcnaao 

— 2nd iteration 

— 3rd iteration 


Fig. 5 Comparison of the testing step response with the iterative LQG 
controller 


suspension system is used to validate this controller design. 
Both numerical simulations and test data show that the control- 
ler converges quickly and is very effective when the system is 
subjected to modeling error and unknown noise statistics. 
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6 Conclusion 

In contrast to most existing LQG controller designs of which 
the great majority solve two separate, but dual problems: the 
LQR and Kalman filter design, this paper proposes an iterative 
LQG controller design approach. A closed-loop identification 
method is developed to update the open-loop state-space model 
and the Kalman filter gain simultaneously from the closed- 
loop input/ output test data. The method is derived under the 
stochastic framework, taking into account the effects of process 
noise as well as measurement noise. For a noise free system, 
the exact open-loop model can be obtained after the first closed- 
loop identification and the identified Kalman filter gain becomes 
the dead-beat observer gain. For a noise corrupted system, itera- 
tions are required to update the open-loop model and the Kal- 
man filter gain from testing until the iterative LQG controller 
converges. In each iteration, since the Kalman filter gain is 
identified directly from test data, the LQG design is simplified 
to be an LQR design. A highly unstable large-angle magnetic 



Fig, 4 Comparison of the (1,1) element of the controller Markov param- 
eters Y d (k) = C d A k a ~'B d 
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APPENDIX 

The analytical model of the large-angle magnetic suspension 
test facility is 

x = A m x + B„u (Al) 

y = C m x (A2) 

where * = { x /),A m = J«], B m = [£*] and C m = [C, 

Osxs]* The state variable x p includes pitch and yaw angles and 
three linear displacements of the cylinder’s centroid. The matri- 
ces A 2 i, A 22 , B 2 , and C, are 


The eigenvalues of the system matrix A m are ±58.78, ±57.81, 
±9.78, ±77.97, and ± 7*0.96. The matrix C t which relates the 
sensor output voltage to the displacement can be obtained from 
calibration and is assumed known. To recover the displacement 
from the sensor output voltage, one can use x p - Cf'y. 

The performance index for the state feedback design is chosen 
as 

oo 

P-/. = £ ylQy k + u T k Ru k (A3) 

*-1 

where Q = (Cr ! ) r diag [l.*3 l.*3 2.*8 2.*8 2.*8]Cr l and R 

= hx s* 


‘3.3415* + 03 0 

0 3.3415* + 03 

-9.8070* + 00 -2.4664* - 15 
-3.6031* - 15 1.9618* - 15 

-2.3357* - 16 -3.6031* - 15 


-3.9392* + 04 
-4.9534* - 12 
4.9937* + 01 
4.3604* - 15 
-2.5089* - 02 


4.9534* - 12 
4.8609* - 12 
4.3604* - 15 
9.5577* + 01 
-9.0007* - 15 


2.0811* - 12 
-1.4472* - 11 
-2.5089* - 02 
-9.0007* - 15 
-9.1324* - 01. 


r 3.8370* + 01 
0 

2.2144* - 01 
0 

-2.7672* - 01 

"8.9024* + 01 
0 

-1.1625* + 02 
0 
0 


3.8370* + 01 
8.9802* + 01 
-1.5274* - 01 
1.2154* - 01 
-8.5465* - 02 


9.5425* + 01 
-1.0725* + 02 


3.8370* + 01 
5.5514* + 01 
7.8453* - 02 
-1.9674* - 01 
2.2388* - 01 

0 

7.8740* + 03 
0 
0 
0 


3.8370* + 01 
-5.5514* + 01 
7.8453* - 02 
1.9674* - 01 
2.2388* - 01 


-6.5359* + 03 
-5.1813* + 03 


3.8370* + 01 
l -8.9802* + 01 
-1.5274* - 01 
-1.2154* - 01 
-8.5465* - 02 

6.0976* + 03" 

0 

6.2500* + 03 . 

0 
0 
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Introduction 

E XTENSIVE research has been conducted in active control 
of dynamic systems. Achieving high control performance on 
these systems usually requires an accurate model. Such a model can 
be derived from system identification techniques using experimen- 
tal data. Recently, a method was developed to identify a state-space 
model from frequency response data for open-loop systems by us- 
ing the state-space frequency domain identification algorithm. 1 This 
method uses a matrix-fraction for the curve fitting, and the curve fit- 
ting is reformulated as a linear problem that can be solved by the 
ordinary least-squares method in one step. 

A different method has been proposed to identify a state-space 
plant from closed-loop I/O time-domain data with known feedback 
dynamics. 2 This Note is an extension of this time-domain closed- 
loop identification method to frequency domain. The method can 
identify a linear open-loop stochastic system from closed-loop fre- 
quency response test data with known feedback dynamics. The 
relationship between the frequency response function (FRF) and 
the closed-loop system and Kalman-filter Markov parameters is de- 
rived for linear stochastic systems. Once the closed-loop system and 
Kalman-filter Markov parameters are obtained from FRF, a recur- 
sive formula for computing the open-loop system and the Kalman- 
filter Markov parameters from the closed-loop system, Kalman filter 
and controller Markov parameters can be used. Finally, the open- 
loop system can be realized from the calculated open-loop system 
Markov parameters. 

Linear State-Space and FRF Relationship 

A finite-dimensional, linear, discrete-time, time-invariant system 
can be modeled as 

•** + i = Ax* -f- Bu k ■+- w k (1) 

y k = Cx k + v k (2) 

where x € R " * 1 , u 6 R 1 * 1 , y e /? m * 1 are state, input, and output 
vectors, respectively; w k is the process noise; v k is the measurement 


noise and [A , B, C] are the state-space parameters. Sequences w k 
and v k are assumed Gaussian, white, and stationary with zero mean 
and covariance matrices Q and R, respectively. One can derive a 
steady-state filter innovation model 1 : 

x k + i = Ax k + Bu k + AKe k (3) 

y k = Cx k + e k (4) 

where x k is the a priori estimated state, K is the steady-state Kalman- 
filter gain, and e k is the residual after filtering: s k = y k - Cx k . On the 
other hand, a dynamic output feedback controller can be modeled 
as 


Pk+ i = A d p k + B d y k (5) 

= Cjp k + D d y k + r k (6) 


where A d , B d , C d , and D d are the system matrices of the controller, 
p 6 R ‘ x 1 is the controller state vector; and r e R * * 1 is the reference 
input to the closed-loop system. Combining Eqs. (3) and (6), the 
augmented closed-loop system dynamics become 


hk+i — A c r] k + B c r k 4- A c K c s k 
yk = C c t] k + e k 
or 

r lk + 1 = Ar} k + B c r k + A c K c y k 

where 

-[:]■ 

C c = [C 0], 


Ac — 


n f D DLj 

B d C A 


* J 


AK + BD d 
B d 


A C K C = 
and A = A c — A c K c C, 


(7) 

( 8 ) 

(9) 


The z transforms of Eqs. (8) and (9) yield 


y(z) - C c t](z) + e(z) (10) 

n(z) = (zl t - Ar'[A c tf c y(z) + B c r{z )] (1 1) 

where I, is an identity matrix with dimension t =n+l. Substituting 
Eq. (11) into Eq. (10), one obtains 

y(z) = [l m — C c (zl, - A)- l A c K c ]~ l C c (zI t - A)~' B c r(z) 

+ [ Im - c c (zl, - A) -1 A c AT c ] -I £(z) ( 12 ) 

The z transforms of the dynamic output feedback controller (5) 
and (6) and the closed-loop state-space model (7) and (8) yield 
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where Y d (k) = C d A k d ~ l B d are the controller Markov parameters, 
Y c (k) = C c A k c ~ l B c are the closed-loop-system Markov parameters, 
and N c (k ) = C c A k ~ 1 A C K C are the closed-loop Kalman-filter Mar- 
kov parameters. Note also that Y d (Q) = D d and N c (0) = I m . 

The transfer-function matrix of the system described by Eqs. (12) 
and (14) is 

G(z~ l )[L - CM - A)- l A e K e y l CM ~ A)~'B C 

PC 

= £ n(«r‘ (is) 

*= 1 

The FRF is simply the transfer function matrix G(z -1 ) calculated 
along the unit circle in the z plane. It is also chosen that the transfer 
function matrix can be expressed by a left-fraction description 1 as 


Table 1 Comparison of eigenvalues of discrete- time 
analytical and identified model 


Analytical model 

Identified 
from simulation 
(1% noise variance) 

Identified 
from testing 

1.1687 

1.1686 

1-2892 

1.1629 

1.1595 

1.2796 

1.0101 

1.0019 

10327 

0.9810 

0.9794 

0.9042 

0.9977 ± 0.0257/ 

0.9051 ±0.0983/ 

1.0094 ± 0.034 1/ 

0.9920 ±0.0133/ 

0.8546 ±0.1947/ 

0-9972 ± 0.0221/ 

0.8633 ± 0.0009/ 

0.8749, 0.9323 

0.8451 ±0.2084/ 


G(z- l )=a~ l (z- l )P(z- 1 ) 


(16) 


where both a(z _l ) and ^(z" 1 ) are matrix polynomials and can be 
found as a solution of the least-squares method. From Eq. (7.19) in 
Ref. 1, one has 


X>-' 


£ ^(Oz" 


i=0 


= £az- 


i=0 



“o = /«. A> = 0 (17) 

From this relation, the closed-loop-system Markov parameters can 
be recursively calculated from the estimated a and /S matrix poly- 
nomials by using the parameter convolution of polynomial products 
as follows: 

k 

K c (*) = ft-£*,T c (*-0 (18) 

i = l 

Similarly, the closed-loop Kalman-filter Markov parameters can be 
recursively calculated from the estimated a matrix polynomials as 
follows: 



Fig. 1 Comparison of the closed-loop analytical (thin line) and recon- 
structed (thick line) input- 1/output-l FRFs. The reconstructed FRF is 
obtained using the identified system matrices. 


k 

N c (k) = ~Yl cti N c(k-i) ( 19 ) 

i = 1 

Then, from the closed-loop-system Markov parameters Y c (k ) and 
the closed-loop Kalman-filter Markov parameters N c (k), one can 
recursively calculate 2 the open-loop-system Markov parameters 
Y(k) — CA k ~ l B and the open-loop Kalman-filter Markov parame- 
ters N (k) ~ CA k ~ 1 A K with the known controller Markov param- 
eters Y d (k) = C d A k ~ l B d , 

j k 

Y(j) = Y e (J) Y ^ Y ^ k ~ i)Y ‘U - k) < 20 > 

k-l i = 1 
j * 

NU) = N c (j) Y ^ Y ^ k - W'U - k) ( 21 > 

*=i / = i 

Note that Y d (0) = D d , N c (0) = / m , and y c (0) = ^ 0 - The open- 
loop state-space model can be realized from the open-loop-system 
Markov parameters through the singular value decomposition 
method. 1 Once the open-loop A and C are obtained, one can 
easily calculate the open-loop Kalman filter gain from the open- 
loop Kalman-filter Markov parameters N(k) in a least-squares 
sense. 2 

Numerical and Test Example 

An example is provided that consists of numerical simulations 
and actual hardware tests to validate the feasibility of the proposed 
frequency-domain closed-loop identification method. The large-gap 
magnetic suspension system 2 consists of a planar array of five cop- 
per electromagnets that actively suspend a small cylinder with a 
permanent magnet core. The cylinder is a rigid body and has six 
independent degrees of freedom, namely, three displacements (x, y, 


and z) and rotations (pitch, yaw, and roll). The roll of the cylinder 
is uncontrollable, and is assumed to be motionless. Because it is 
difficult to accurately model the magnetic field and its gradients, the 
analytical model needs to be improved through identification from 
experimental data. The discrete-time state-space parameters of the 
system and the dynamic output feedback using a sampling rate of 
250 Hz are shown in the Appendix. 

Table I compares identified system eigenvalues with true ones 
from numerical simulations. The results show perfect match when 
there is no noise and quite good agreement even with 1% of pro- 
cessing and measurement noises (1% noise variance). Figure 1 also 
shows the comparison of the closed-loop analytical (thin line) and 
reconstructed (thick line) input- 1/output-l FRFs v 

Five experiments also were performed. In each experiment, only 
one of the five actuators had a single random reference input, and 
all others had zero reference input. A total of 4096 data points at 
a sampling rate of 250 Hz from each sensor were recorded. Six 
FRFs from these single input/six output data can be derived. The 
experiment is repeated by sending the same single random input to a 
different actuator each time. Thirty FRFs were obtained. The order 
of the matrix polynomial was set to 13. The identified eigenvalues 
from testing are shown in Table 1. 

Concluding Remarks 

A method of identifying a linear state-space model of a plant 
from closed-loop frequency response data with known feedback 
dynamics is developed. The main contribution is the derivation of 
the relationship between the open-loop-system Markov parame- 
ters and the closed-loop frequency response function for a linear 
stochastic system. Numerical simulations and experimental results 
of a highly unstable large-gap magnetic suspension system are pre- 
sented. 
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Appendix: Analytical System and Controller Matrices 

The matrices of Eqs. (1), (2), (5), and (6) are as follows: 


' 1.1687 
-0.0000 
-0.0000 
-0.0000 
0.0000 
0.0000 
0.0000 
- 0.0000 
0.0000 
0.0000 


0.0006 

1.1629 

0.0001 

0.0000 

0.0002 

-0.0000 

- 0.0001 

- 0.0000 

- 0.0000 

-0.0000 


- 0.0000 

-0.0000 

1.0178 

0.0001 

-0.0004 

- 0.0021 

-0.0064 

0.0109 

-0.0086 

0.0004 


0.0000 

-0.0000 

-0.0017 

1.0051 

0.0008 

-0.0240 

- 0.0001 

-0.0009 

0.0009 

0.0002 


0.0000 

-0.0000 

-0.0037 

0.0001 

1.0106 

0.0005 

-0.0213 

-0.0045 

0.0032 

0.0006 


0.0000 

-0.0000 

0.0021 

0.0295 

-0.0018 

0.9908 

-0.0041 

0.0021 

0.0009 

0.0012 


0.0000 

0.0000 

0.0074 

0.0006 

0.0223 

0.0028 

0.9692 

0.0050 

0.0031 

0.0545 


-0.0000 

-0.0000 

-0.0127 

0.0015 

0.0066 

- 0.0010 

0.0064 

0.9260 

-0.0589 

- 0.0002 


-0.0000 

0.0000 

0.0112 

- 0.0011 

-0.0039 

0.0003 

0.0004 

-0.0549 

0.9125 

- 0.0002 


-0.0000' 
-0.0000 
0.0006 
0.0003 
0.0030 
- 0.0011 
0.0003 
0.0028 
-0.0008 
0.8652 I 


" 0.0035 

0.0706 

0.0519 

-0.0363 

-0.0434 

-0.0326 

-0.0340 

-0.0425 

0.0580 

-0.0454 

0.0983 

-0.0361 

-0.0926 

-0.0315 

0.0881 

0.0865 

0.1160 

0.0124 

0.0263 

0.0982 

-0.1015 

-0.0368 

0.1033 

0.0854 

0.1373 

0.0057 

0.0719 

0.0859 

-0.0159 

-0.0637 

-0.1326 

0.1165 

0.0158 

-0.1531 

-0.0261 

0.0041 

_ -0.0484 

-0.0800 

-0.0513 

-0.0553 


‘ -0.03 13 0.4029 -0.0469 0.2269 -0.0381 

0.0291 -0.4213 0.0006 0.2248 0.0290 

-0.4423 0.1071 0.1809 0.0553 0.0669 

-0.4254 -0.1184 -0.1787 -0.0092 -0.0829 

0.4495 -0.0763 0.0574 0.0273 -0.1861 

0.3889 0.1015 -0.0614 0.0085 0.1739 


-0.4423 0.1071 


-0.0633 

-0.0396 

0.0254 

-0.0218 

-0.0242 

-0.0154 

-0.0066 

0.0625 

0.1245 

-0.1009 


0.2269 -0.0381 -0.1961 0.1274 -0.0363 


0.0290 -0.2097 -0.1079 -0.0130 


0.0198 

0.0297 


-0.1513 

0.1502 


0.0669 -0.0618 -0.0906 -0.0418 -0.2228 -0.0472 


0.0200 


0.0273 -0.1861 -0.0400 0.1239 

0.0085 0.1739 0.0012 -0.1277 


0.1217 -0.2197 -0.0559 
0.1239 0.2109 0.0827 

0.1277 0.0386 0.1913 


■0.0559 0.0630 

0.0827 0.0464 

0.1913 -0.0634 


0.3333 

0 

0 

0 

0 


0 

0.3333 


0 0.6000 


0 0.6000 

0 0 


-0.0206 0.0206 
0 0 

0 0 

0 0 


0 0.6000 J 


L 0.0004 

0.0004 


"0.0796 

0.0000 

7.3872 

0.0000 


0.1032 

0.0716 

-5.9772 

4.3222 

.0e + 03* 

0.0886 

0.0442 

2.2836 

-6.9917 


0.0886 

-0.0442 

2.2836 

6.9917 


_ 0.1 032 

-0.0716 

-5.9772 

-4.3222 


10.8171 

3.9903 

-7.0133 

7.0133 

7.0133 

6.7151 

-2.1362 

11.2687 

-8.3349 

-3.0144 

-2.1923 

-9.7904 

-7.9381 

9.7505 

-5.4144 

-2.1923 

-9.7904 

3.6020 

-5.4144 

9.7505 

6.7151 

-2.1362 

0.0807 

-3.0144 

-8.3349 


0 0 0 0 
-0.0098 -0.0098 0.0098 0.0098 

0.0003 -0.0003 -0.0003 0.0003 

-0.0003 0.0003 -0.0003 0.0003 

0 0 0 0 


-5.5493 ‘ 
-1.7160 
4.4907 
4.4907 
— 1.7160 _ 

— 7.0133 
0.0807 
3.6020 
-7.9381 
11.2687 
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